Preparación de datos

library(tidyverse)

## ── Attaching packages ─────────────────────────────────────── tidyverse 1.3.1 ──

## ✓ ggplot2 3.3.3     ✓ purrr   0.3.4
## ✓ tibble  3.1.1     ✓ dplyr   1.0.5
## ✓ tidyr   1.1.3     ✓ stringr 1.4.0
## ✓ readr   1.4.0     ✓ forcats 0.5.1

## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## x dplyr::filter() masks stats::filter()
## x dplyr::lag()    masks stats::lag()

library(tidymodels)

## Registered S3 method overwritten by 'tune':
##   method                   from   
##   required_pkgs.model_spec parsnip

## ── Attaching packages ────────────────────────────────────── tidymodels 0.1.3 ──

## ✓ broom        0.7.6      ✓ rsample      0.1.0 
## ✓ dials        0.0.10     ✓ tune         0.1.6 
## ✓ infer        1.0.0      ✓ workflows    0.2.3 
## ✓ modeldata    0.1.1      ✓ workflowsets 0.1.0 
## ✓ parsnip      0.1.7      ✓ yardstick    0.0.8 
## ✓ recipes      0.1.16

## ── Conflicts ───────────────────────────────────────── tidymodels_conflicts() ──
## x scales::discard() masks purrr::discard()
## x dplyr::filter()   masks stats::filter()
## x recipes::fixed()  masks stringr::fixed()
## x dplyr::lag()      masks stats::lag()
## x yardstick::spec() masks readr::spec()
## x recipes::step()   masks stats::step()
## ● Use tidymodels_prefer() to resolve common conflicts.

library(knitr)
library(kableExtra)

## 
## Attaching package: 'kableExtra'

## The following object is masked from 'package:dplyr':
## 
##     group_rows

library(corrr)
library(GGally)

## Registered S3 method overwritten by 'GGally':
##   method from   
##   +.gg   ggplot2

library(ggplot2)

dt_train = read_csv("encuesta_salud_train.csv")

## 
## ── Column specification ────────────────────────────────────────────────────────
## cols(
##   record = col_double(),
##   edad = col_double(),
##   genero = col_character(),
##   nivel_educativo = col_character(),
##   altura = col_double(),
##   peso = col_double(),
##   frecuencia_hambre_mensual = col_character(),
##   dias_consumo_comida_rapida = col_double(),
##   edad_consumo_alcohol = col_character(),
##   consumo_diario_alcohol = col_double(),
##   dias_actividad_fisica_semanal = col_double(),
##   consumo_semanal_frutas = col_character(),
##   consumo_semanal_verdura = col_character(),
##   consumo_semanal_gaseosas = col_character(),
##   consumo_semanal_snacks = col_character(),
##   consumo_semanal_comida_grasa = col_character()
## )

dt_test = read_csv("encuesta_salud_test.csv")

## 
## ── Column specification ────────────────────────────────────────────────────────
## cols(
##   record = col_double(),
##   edad = col_double(),
##   genero = col_character(),
##   nivel_educativo = col_character(),
##   altura = col_double(),
##   peso = col_double(),
##   frecuencia_hambre_mensual = col_character(),
##   dias_consumo_comida_rapida = col_double(),
##   edad_consumo_alcohol = col_character(),
##   consumo_diario_alcohol = col_double(),
##   dias_actividad_fisica_semanal = col_double(),
##   consumo_semanal_frutas = col_character(),
##   consumo_semanal_verdura = col_character(),
##   consumo_semanal_gaseosas = col_character(),
##   consumo_semanal_snacks = col_character(),
##   consumo_semanal_comida_grasa = col_character()
## )

dt_train_6 = read_csv("encuesta_salud_modelo6.csv")

## 
## ── Column specification ────────────────────────────────────────────────────────
## cols(
##   record = col_double(),
##   edad = col_double(),
##   genero = col_character(),
##   nivel_educativo = col_character(),
##   altura = col_double(),
##   peso = col_double(),
##   frecuencia_hambre_mensual = col_character(),
##   dias_consumo_comida_rapida = col_double(),
##   edad_consumo_alcohol = col_character(),
##   consumo_diario_alcohol = col_double(),
##   dias_actividad_fisica_semanal = col_double(),
##   consumo_semanal_frutas = col_character(),
##   consumo_semanal_verdura = col_character(),
##   consumo_semanal_gaseosas = col_character(),
##   consumo_semanal_snacks = col_character(),
##   consumo_semanal_comida_grasa = col_character()
## )

#Se renombran variables en los datasets de training y testing para facilitar la lectura-
dt_train <- dt_train %>% rename( nivel_edu = nivel_educativo,
                     frec_hambre_m = frecuencia_hambre_mensual,
                     q_comida_rapida = dias_consumo_comida_rapida,
                     edad_alcohol = edad_consumo_alcohol,
                     frec_alcohol_d =consumo_diario_alcohol,
                     act_fisica_s = dias_actividad_fisica_semanal,
                     frec_frutas_s = consumo_semanal_frutas,
                     frec_verdura_s = consumo_semanal_verdura,
                     frec_gaseosas_s = consumo_semanal_gaseosas,
                     frec_snacks_s = consumo_semanal_snacks,
                     frec_comida_grasa_s = consumo_semanal_comida_grasa
                     )

dt_test <- dt_test %>% rename( nivel_edu = nivel_educativo,
                     frec_hambre_m = frecuencia_hambre_mensual,
                     q_comida_rapida = dias_consumo_comida_rapida,
                     edad_alcohol = edad_consumo_alcohol,
                     frec_alcohol_d =consumo_diario_alcohol,
                     act_fisica_s = dias_actividad_fisica_semanal,
                     frec_frutas_s = consumo_semanal_frutas,
                     frec_verdura_s = consumo_semanal_verdura,
                     frec_gaseosas_s = consumo_semanal_gaseosas,
                     frec_snacks_s = consumo_semanal_snacks,
                     frec_comida_grasa_s = consumo_semanal_comida_grasa
                     )

dt_train_6 <- dt_train_6 %>% rename( nivel_edu = nivel_educativo,
                     frec_hambre_m = frecuencia_hambre_mensual,
                     q_comida_rapida = dias_consumo_comida_rapida,
                     edad_alcohol = edad_consumo_alcohol,
                     frec_alcohol_d =consumo_diario_alcohol,
                     act_fisica_s = dias_actividad_fisica_semanal,
                     frec_frutas_s = consumo_semanal_frutas,
                     frec_verdura_s = consumo_semanal_verdura,
                     frec_gaseosas_s = consumo_semanal_gaseosas,
                     frec_snacks_s = consumo_semanal_snacks,
                     frec_comida_grasa_s = consumo_semanal_comida_grasa
                     )

Análisis exploratorio

Estructura y variables

El archivo de training tiene 7024 filas y 16 columnas. Las variables con las que se van a trabajar contienen los siguientes datos:

consumo_diario_alcohol (medido en tragos): (0.0, 0.5, 1.0, 2.0, 3.0, 4.0 y 5.0).
dias_actividad_fisica_semanal: (0, 1, 2, 3, 4, 5, 6, 7)
dias_consumo_comida_rapida: (0, 1, 2, 3, 4, 5, 6, 7)
edad: (12, 13, 14, 15, 16, 17, 18)
altura: (de 125 a 200)
peso: (de 27 a 144)
genero: (Femenino“,”Masculino)
nivel_educativo: (1er año/10mo grado nivel Polimodal o 3er año nivel Secundario, 2do año/11vo grado nivel Polimodal o 4to año nivel Secundario, 3er año/12vo grado nivel Polimodal o 5to año nivel Secundario, 8vo grado nivel Primario/Polimodal o 1er año nivel Secundario,9no grado nivel Primario/Polimodal o 2do año nivel Secundario, Dato perdido).
frecuencia_hambre_mensual: (Algunas veces, Casi siempre, Dato perdido, Nunca, Rara vez, Siempre)
edad_consumo_alcohol: (10 o 11 años, 12 o 13 años, 14 o 15 años, 16 o 17 años, 18 años o más, 7 años o menos, 8 o 9 años, Nunca tomé alcohol más que unos pocos sorbos)
consumo_semanal_frutas: (1 a 3 veces durante los últimos 7 días, 1 vez al día, 2 veces al día, 3 veces al día, 4 a 6 veces durante los últimos 7 días, 4 o más veces al día, Dato perdido, No comí frutas durante los últimos 7 días)
consumo_semanal_verdura: (1 a 3 veces durante los últimos 7 días, 1 vez al día, 2 veces al día, 3 veces al día, 4 a 6 veces durante los últimos 7 días, 4 o más veces al día, Dato perdido, No comí verduras ni hortalizas durante los últimos 7 días)
consumo_semanal_gaseosas: (1 a 3 veces durante los últimos 7 días, 1 vez al día, 2 veces al día, 3 veces al día, 4 a 6 veces durante los últimos 7 días, 4 o más veces al día, Dato perdido, No tomé gaseosas en los últimos 7 días)
consumo_semanal_snacks: (1 a 3 veces durante los últimos 7 días, 1 vez al día, 2 veces al día, 3 veces al día, 4 a 6 veces durante los últimos 7 días, 4 o más veces al día, Dato perdido, No comí comida salada o snacks en los últimos 7 días)
consumo_semanal_comida_grasa: (1 a 3 veces durante los últimos 7 días, 1 vez al día, 2 veces al día, 3 veces al día, 4 a 6 veces durante los últimos 7 días, 4 o más veces al día, Dato perdido, No comí comida alta en grasa en los últimos 7 días)

Como variables numéricas están: altura y peso. Luego consumo_diario_alcohol, dias_actividad_fisica_semanal, dias_consumo_comida_rapida y edad contienen valores discretos y acotados al dominio, por lo que se podrían considerar categóricas. El resto de las variables son todas categóricas.

dt_train %>% glimpse()

## Rows: 7,024
## Columns: 16
## $ record              <dbl> 502, 26488, 31473, 14154, 36578, 53730, 11892, 355…
## $ edad                <dbl> 17, 15, 15, 16, 17, 15, 13, 17, 17, 16, 16, 14, 15…
## $ genero              <chr> "Femenino", "Masculino", "Masculino", "Masculino",…
## $ nivel_edu           <chr> "2do año/11vo grado nivel Polimodal o 4to año nive…
## $ altura              <dbl> 165, 178, 172, 170, 170, 178, 156, 163, 164, 167, …
## $ peso                <dbl> 62, 62, 62, 65, 75, 88, 46, 60, 57, 51, 100, 33, 6…
## $ frec_hambre_m       <chr> "Rara vez", "Rara vez", "Nunca", "Nunca", "Rara ve…
## $ q_comida_rapida     <dbl> 0, 0, 3, 1, 1, 2, 0, 0, 0, 3, 4, 2, 1, 1, 3, 0, 0,…
## $ edad_alcohol        <chr> "14 o 15 años", "7 años o menos", "Nunca tomé alco…
## $ frec_alcohol_d      <dbl> 5.0, 4.0, 0.0, 0.0, 0.0, 5.0, 1.0, 0.5, 5.0, 0.0, …
## $ act_fisica_s        <dbl> 7, 7, 7, 7, 0, 7, 0, 2, 7, 3, 2, 2, 7, 1, 4, 0, 1,…
## $ frec_frutas_s       <chr> "No comí frutas durante los últimos 7 días", "No c…
## $ frec_verdura_s      <chr> "4 a 6 veces durante los últimos 7 días", "4 a 6 v…
## $ frec_gaseosas_s     <chr> "1 a 3 veces durante los últimos 7 días", "1 a 3 v…
## $ frec_snacks_s       <chr> "1 a 3 veces durante los últimos 7 días", "No comí…
## $ frec_comida_grasa_s <chr> "No comí comida alta en grasa en los últimos 7 día…

Correlación entre variables numéricas

dt_train %>% select(edad, altura, peso,frec_alcohol_d, act_fisica_s, q_comida_rapida) %>%
  correlate() %>% rplot()

## 
## Correlation method: 'pearson'
## Missing treated using: 'pairwise.complete.obs'

## Don't know how to automatically pick scale for object of type noquote. Defaulting to continuous.

# Correlación entre variables numéricas (las más significativas)
cor(dt_train$peso, dt_train$edad)

## [1] 0.2846795

cor(dt_train$peso, dt_train$altura)

## [1] 0.5759839

cor(dt_train$peso, dt_train$q_comida_rapida)

## [1] -0.05872582

cor(dt_train$peso, dt_train$frec_alcohol_d)

## [1] 0.07828881

cor(dt_train$peso, dt_train$act_fisica_s)

## [1] 0.06693092

cor(dt_train$altura, dt_train$edad)

## [1] 0.2539131

cor(dt_train$frec_alcohol_d, dt_train$edad)

## [1] 0.2288491

En cuanto a las correlaciones, la más alta es entre peso y altura con 0.5759839, luego le sigue la edad y peso con 0.2846795. Finalmente altura y edad con 0.2539131 y frec_alcohol_d y edad con 0.2288491, siendo todas asociaciones lineales positivas.

Apertura por género

dt_train %>% 
  select(genero,edad, altura, frec_alcohol_d, act_fisica_s, q_comida_rapida, peso) %>% 
  mutate(genero = factor(genero)) %>%  
  ggpairs(., mapping = aes(colour = genero), title = "Matriz de correlaciones",
          upper = list(continuous = wrap("cor", size = 3, hjust=0.6)), legend = 25, progress=FALSE) + 
  theme(axis.text.x = element_text(angle=45, vjust=0.5), legend.position = "bottom") +
  theme_bw()

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

ggplot(data = dt_train, aes(x = peso, y= altura, color = genero)) +
  geom_point(alpha = 0.75) + 
  labs(title = "Peso y Altura por género") + 
  theme(legend.position = 'none') + 
  facet_wrap(~ genero)

ggplot(data = dt_train, aes(x = peso, y= edad, color = genero)) +
  geom_point(alpha = 0.75) + 
  labs(title = "Peso y Edad por género") + 
  theme(legend.position = 'none') + 
  facet_wrap(~ genero)

ggplot(data = dt_train, aes(x = altura, y= edad, color = genero)) +
  geom_point(alpha = 0.75) + 
  labs(title = "Altura y Edad por género") + 
  theme(legend.position = 'none') + 
  facet_wrap(~ genero)

ggplot(data = dt_train, aes(x = frec_alcohol_d, y= edad, color = genero)) +
  geom_point(alpha = 0.75) + 
  labs(title = "Frecuencia de alcohol diario y Edad por género") + 
  theme(legend.position = 'none') + 
  facet_wrap(~ genero)

En el caso de la relación entre peso y altura según el género, se observa que los hombres alcanzan mayores alturas y mayores pesos que las mujeres. Además en algunos casos las mujeres tienen mayor peso que los hombres dada la misma altura. Luego, como se comentó anteriormente, tanto la edad como la frecuencia_alcohol_diaria están discretizadas, viendose barras más que una nube de puntos.

Hambre mensual, consumo semanal de verdura y consumo semanal de comida grasa

ggplot(dt_train, aes(x= frec_hambre_m, fill= frec_verdura_s)) +
  geom_bar(aes(y = (..count..)/sum(..count..))) +
  labs(title = "Distribución de Frecuencia de Hambre",
       subtitle = "(Subdividido por consumo de verdura)") + 
  guides(fill = guide_legend(title = "Verduras")) +
  scale_x_discrete("Hambre") + 
  scale_y_continuous("", labels=scales::percent) +
  theme(axis.text.x = element_text(angle = 70, hjust = 1))

ggplot(dt_train, aes(x= frec_hambre_m, fill= frec_comida_grasa_s)) +
  geom_bar(aes(y = (..count..)/sum(..count..))) +
  labs(title = "Distribución de Frecuencia de Hambre",
       subtitle = "(Subdividido por consumo de comida grasa)") + 
  guides(fill = guide_legend(title = "Comidas grasas")) +
  scale_x_discrete("Hambre") + 
  scale_y_continuous("", labels=scales::percent) +
  theme(axis.text.x = element_text(angle = 70, hjust = 1))

La gran mayoría, “nunca” paso hambre, siguiéndoles los que “rara vez” pasaron hambre y finalmente los que “algunas veces” pasaron hambre. Las demás categorías (“casi siempre”, “dato perdido” y “siempre”) contienen muy pocas personas.

Ahora bien, si observamos dichas categorías relacionadas al consumo de verduras, se observa que gran parte consume verduras en mayor o menor medida, siendo la opción ganadora de “1 a 3 veces en los últimos 7 días” y siguiéndole de “4 a 6 veces en los últimos 7 días” para las categorías “nunca” y “rara vez”, mientras que para “algunas veces” hay un empate entre ellas.

Por último, si observamos dichas categorías relacionadas al consumo de comidas grasas, vemos que para las categorías “nunca”, “rara vez” y “algunas veces” la opción ganadora de consumo de grasas se da de “1 a 3 veces en los últimos 7 días”. Luego en el caso de “nunca” la opción siguiente es “no comió grasas en los últimos 7 días”, mientras que en “rara vez” y “algunas veces” hay un empate entre “4 a 6 veces en los últimos 7 días” versus “no comió grasas en los últimos 7 días”.

Finalmente, gran parte de las personas comen igual frecuencia de verduras y comidas grasas (1 a 3 veces durante los últimos 7 días).

Modelo Inicial

modelo_inicial = lm(formula = peso ~ altura + edad + genero + act_fisica_s + frec_alcohol_d, data= dt_train)
tidy_modelo_inicial <- tidy(modelo_inicial, conf.int = TRUE)
tidy_modelo_inicial

Coeficientes estimados

Interpretación

Intercept: (Categoría basal de la variable categórica) El valor del coeficiente estimado es la media del peso para personas sin altura, sin edad con 0 actividad física semanal y con 0 consumo de alcohol diario. En este modelo es -68.922688070 kilogramos, lo cual no tiene sentido ya que una persona debería tener al menos altura y edad. Además no podría tener un peso negativo.
generoMasculino: El valor del coeficiente estimado es la diferencia en los niveles medios del género masculino respecto al género femenino (categoría basal). Es decir, este coeficiente = 1.262643558 indica cuanto más alto es el peso para el género masculino respecto del género femenino (categoría basal), dado una altura, una edad,una actividad física semanal y un consumo de alcohol diario.
Altura: El coeficiente estimado es 0.650606544 kilogramos, lo que indica que si mantenemos la edad, el género,la actividad física semanal y el consumo de alcohol diario constantes, cada incremento adicional de un centímetro de altura corresponde a un aumento de 0.650606544 kilogramos, en promedio en el peso de la persona.
Edad: El coeficiente estimado es 1.406727060 kilogramos, lo que indica que si mantenemos la altura, el genero, la actividad física semanal y el consumo de alcohol diario constantes, cada incremento adicional de un año en la edad corresponde a un aumento de 1.406727060 kilogramos, en promedio en el peso de la persona.
act_fisica_s: El coeficiente estimado es -0.087391031 kilogramos, lo que indica que si mantenemos la altura, el genero,la edad y el consumo de alcohol diario constantes, cada incremento adicional de un día de actividad física corresponde a una disminución de 0.087391031 kilogramos, en promedio en el peso de la persona.
frec_alcohol_d: El coeficiente estimado es 0.007271379 kilogramos, lo que indica que si mantenemos la altura, el genero, la edad y la actividad física semanal constantes, cada incremento adicional de un trago corresponde a un aumento de 0.007271379 kilogramos, en promedio en el peso de la persona.

Significatividad

Se observa que las variables altura, edad y genero masculino resultan estadísticamente significativas para explicar el peso de las personas (p-valores < 0.05, resultado de aplicar el test t). Además, se puede apreciar que los intervalos de confianza (IC) del 95% de dichas variables no contienen al 0. Por otro lado, las variables actividad física semanal y consumo de alcohol diario no resultan estadísticamente significativas (p-valores > 0.05 y sus ICs conteniendo al 0).

La tabla de ANOVA muestra que, según el resultado del test F, la variable genero en su conjunto resulta estadísticamente significativa para explicar el peso (p-valor < 0.05)

tidy(anova(modelo_inicial))

Observando el resultado del test-F (Test de significatividad global) vemos que el modelo si resulta significativo para explicar el peso (P-valor < 0.05).
Además, el modelo logra explicar una variablidad del 35.2 %, observando el R² ajustado

glance(modelo_inicial)

Modelo categóricas

dt_train_1 = dt_train
#Se cambia la categoría basal
dt_train_1 <- dt_train_1 %>% mutate(frec_snacks_s = as.factor(frec_snacks_s))
dt_train_1 <- dt_train_1 %>% mutate(frec_snacks_s = relevel(frec_snacks_s, ref = "No comí comida salada o snacks en los últimos 7 días"))

#Se crea un modelo con interación entre edad y género.
modelo_interacion = lm(formula = peso ~ altura + frec_snacks_s + genero*edad, data= dt_train_1)
tidy_modelo_interaccion <- tidy(modelo_interacion, conf.int = TRUE) %>%  arrange(p.value)
tidy_modelo_interaccion

Coeficientes estimados

Interpretación

consumoSemanalSnacks

frec_snacks_s1 a 3 veces durante los últimos 7 días: El coeficiente estimado muestra cuando disminuye (1.35163 kilogramos) el peso medio para aquellas personas que comen snacks de 1 a 3 veces respecto de aquellos que no comieron snacks en los ultimos 7 días (categoría basal), dada la altura, edad y genero.
frec_snacks_s1 vez al día: El coeficiente estimado muestra cuando disminuye (0.60788 kilogramos) el peso medio para aquellas personas que comen snacks de 1 vez al día respecto de aquellos que no comieron snacks en los ultimos 7 días (categoría basal), dada la altura, edad y genero.
frec_snacks_s2 veces al día: El coeficiente estimado muestra cuando disminuye (1.09466 kilogramos) el peso medio para aquellas personas que comen snacks de 2 veces al día respecto de aquellos que no comieron snacks en los ultimos 7 días (categoría basal), dada la altura, edad y genero.
frec_snacks_s3 veces al día: El coeficiente estimado muestra cuando disminuye (1.27596 kilogramos) el peso medio para aquellas personas que comen snacks de 3 veces al día respecto de aquellos que no comieron snacks en los ultimos 7 días (categoría basal), dada la altura, edad y genero.
frec_snacks_s4 a 6 veces durante los últimos 7 días: El coeficiente estimado muestra cuando disminuye (2.27004 kilogramos) el peso medio para aquellas personas que comen snacks de 4 a 6 veces durante los últimos 7 días respecto de aquellos que no comieron snacks en los ultimos 7 días (categoría basal), dada la altura, edad y genero.
frec_snacks_s4 o más veces al día: : El coeficiente estimado muestra cuando disminuye (2.56697 kilogramos) el peso medio para aquellas personas que comen snacks de 4 o más veces al día respecto de aquellos que no comieron snacks en los ultimos 7 días (categoría basal), dada la altura, edad y genero.
frec_snacks_sDato perdido: El coeficiente estimado muestra cuando disminuye (4.44042 kilogramos) el peso medio para aquellas personas que no completaron la opción respecto de aquellos que no comieron snacks en los ultimos 7 días (categoría basal), dada la altura, edad y genero.

En el caso de la variable consumoSemanalSnacks, se observa que algunos tipos de consumo de snacks resultan estadísticamente significativos para explicar el peso (p-valores < 0.05) mientras otros tipos de consumo no ( frec_snacks_s1 vez al día, frec_snacks_s2 veces al día y frec_snacks_s3 veces al día). Lo mismo se observa a través de los intervalos de confianza del 95%, donde algunos contienen al cero y otros no. Se chequea si los valores medios del peso de las personas son los mismos para los distintos tipos de consumo de snacks respecto del tipo “No comí comida salada o snacks en los últimos 7 días” (categoría basal).

genero.edad

Como la variable generoMasculino asume solamente valores 0 y 1, el termino de la interacción generoMasculino:edad valdrá 0 si generoMasculino = 0 (o sea para las mujeres) y será igual a edad siempre que generoMasculino = 1 (o sea para los hombres).

edad: El coeficiente estimado es la pendiente de la edad en el grupo de mujeres. Indica que por cada aumento en 1 unidad de la edad entre las mujeres, el peso medio aumenta 1.22 unidades.
generoMasculino:edad: El coeficiente estimado representa el aumento en 0.39 de la pendiente en el grupo de hombres con respecto al de las mujeres.

La interacción generoMasculino:edad resulta estadísticamente significativa (p-valor < 0.05), es decir que la edad tiene un efecto diferente en el peso dependiendo del género de la persona. Notar que la variable generoMasculino deja de ser significativa.

Este modelo logra explicar una variablidad del 35.6%

glance(modelo_interacion)

Variable consumoSemanalSnacks

Evaluación de significatividad conjunta

tidy(anova(modelo_interacion))

La tabla de ANOVA muestra que, según el resultado del test F, la variable frec_snacks_s en su conjunto resulta estadísticamente significativa para explicar al peso (p-valor < 0.05). Es decir, que pese a que algunas categorías en su comparación individual con la categoría basal sean poco significativas, la variable en su conjunto sí resulta significativa para el modelo.

Redefinición de las categorías:

Teniendo en cuenta que las siguientes categorías no son significativas: frec_snacks_s1 vez al día, frec_snacks_s2 veces al día y frec_snacks_s3 veces al día, se propone un nuevo agrupamiento para la variable frec_snacks_s. Luego, todas las categorías siguen igual, salvo frec_snacks_s1 vez al día, frec_snacks_s2 veces al día, frec_snacks_s3 veces al día y 4 a 6 veces durante los últimos 7 días que son incluidas en una única nueva categoría “Todos los días en los últimos 7 días”. La clasificación original usaba por un lado como referencia los “ultimos 7 días”, pero además incluía demasiado detalle al suponer que la persona comía todos los días snacks (1, 2, 3, 4 o 6 veces al día durante los ultimos 7 días). Entonces se decide una categorización más homogenea.

nueva_categoria = "Todos los días en los últimos 7 días"
# Se agrega la nueva agrupación de categorías en una nueva variable "frec_snacks_s_v2", tanto para el dataset de training como para el dataset de testing.

dt_train_2 <- dt_train_1 %>% 
  mutate(frec_snacks_s_v2 = case_when(frec_snacks_s == "1 vez al día" ~ nueva_categoria,
                frec_snacks_s == "2 veces al día" ~ nueva_categoria,
                frec_snacks_s == "3 veces al día" ~ nueva_categoria,
                frec_snacks_s == "4 o más veces al día" ~ nueva_categoria,
                frec_snacks_s == "No comí comida salada o snacks en los últimos 7 días" ~  "No comí comida salada o snacks en los últimos 7 días",
                frec_snacks_s == "1 a 3 veces durante los últimos 7 días" ~ "1 a 3 veces durante los últimos 7 días",
                frec_snacks_s == "4 a 6 veces durante los últimos 7 días" ~ "4 a 6 veces durante los últimos 7 días" ,
                frec_snacks_s == "Dato perdido" ~ "Dato perdido" ))


dt_test_2 <- dt_test %>% 
  mutate(frec_snacks_s_v2 = case_when(frec_snacks_s == "1 vez al día" ~ nueva_categoria,
                frec_snacks_s == "2 veces al día" ~ nueva_categoria,
                frec_snacks_s == "3 veces al día" ~ nueva_categoria,
                frec_snacks_s == "4 o más veces al día" ~ nueva_categoria,
                frec_snacks_s == "No comí comida salada o snacks en los últimos 7 días" ~  "No comí comida salada o snacks en los últimos 7 días",
                frec_snacks_s == "1 a 3 veces durante los últimos 7 días" ~ "1 a 3 veces durante los últimos 7 días",
                frec_snacks_s == "4 a 6 veces durante los últimos 7 días" ~ "4 a 6 veces durante los últimos 7 días" ,
                frec_snacks_s == "Dato perdido" ~ "Dato perdido" ))

Notar que ahora cada categoria de frec_snacks_s2 resulta estadisticamente significativa (p-valores < 0.05)

#Se cambia la categoria basal
dt_train_2 <- dt_train_2 %>% 
  mutate(frec_snacks_s_v2 = as.factor(frec_snacks_s_v2))  %>% 
  mutate(frec_snacks_s_v2 = relevel(frec_snacks_s_v2, ref = "No comí comida salada o snacks en los últimos 7 días"))

modelo_interacion_2 = lm(formula = peso ~ altura + frec_snacks_s_v2 + genero*edad, data= dt_train_2)
tidy_modelo_interaccion_2 <- tidy(modelo_interacion_2, conf.int = TRUE) %>%  arrange(p.value)
tidy_modelo_interaccion_2

Luego, se mantiene la variabilidad explicada en 35.6%.

glance(modelo_interacion_2)

Modelos propios y evaluación

Creación de 2 modelos múltiples adicionales

Modelo 1 propuesto

Se decide probar con todas las variables para aprovechar toda la información con que se cuenta.

modelo_propuesto_1 = lm(formula = peso ~ altura + edad + genero + nivel_edu + frec_hambre_m + q_comida_rapida + edad_alcohol + frec_alcohol_d + act_fisica_s + frec_frutas_s + frec_verdura_s + frec_gaseosas_s + frec_snacks_s + frec_comida_grasa_s, data= dt_train)

Modelo 2 propuesto

Se observa que las personas que nunca pasaron hambre son 4844 mientras que el resto de las categorías suma 2180, por eso se decide crear la variable paso_hambre indicando si paso o no hambre en el último mes.

dt_train %>% group_by(frec_hambre_m) %>% count()

#Cambiar la variable "frec_hambre_m" por "paso_hambre" (pasó hambre en el último mes) .
dt_train_propuesta_2 = dt_train %>% 
  mutate(paso_hambre = case_when(frec_hambre_m == "Nunca" ~ "No",
                                 frec_hambre_m == "Algunas veces" ~ "Si",
                                 frec_hambre_m == "Casi siempre" ~ "Si",
                                 frec_hambre_m == "Dato perdido" ~ "Dato perdido",
                                 frec_hambre_m == "Rara vez" ~ "Si",
                                 frec_hambre_m == "Siempre" ~ "Si"))

dt_test_propuesta_2 = dt_test %>% 
  mutate(paso_hambre = case_when(frec_hambre_m == "Nunca" ~ "No",
                                 frec_hambre_m == "Algunas veces" ~ "Si",
                                 frec_hambre_m == "Casi siempre" ~ "Si",
                                 frec_hambre_m == "Dato perdido" ~ "Dato perdido",
                                 frec_hambre_m == "Rara vez" ~ "Si",
                                 frec_hambre_m == "Siempre" ~ "Si"))

modelo_propuesto_2 = lm(formula = peso ~ altura + edad + genero + paso_hambre, data= dt_train_propuesta_2)

Evaluación de la performance de los modelos

En training

#Se arman listas de los modelos
modelos <- list(modelo_ini = modelo_inicial, modelo_interac_2 = modelo_interacion_2, modelo_prop_1 = modelo_propuesto_1, modelo_prop_2 = modelo_propuesto_2)
#Se calculan las métricas para todos los modelos
purrr::map_df(modelos, broom::glance, .id = "modelo") %>% arrange(desc(adj.r.squared))

De los 4 modelos, el que mejor permite medir el porcentaje de variabilidad explicada del peso es el modelo propuesto 1 con un R² ajustado de 0.3608120 (aunque su diferencia con el siguiente modelo categóricas es solo de 0.0033).

#Se aplica función augment -para predecir el peso- para los modelos incluidos en la lista.
lista_pred_training = map(.x = modelos, .f = augment)
#Se obtiene el R cuadrado, RMSE y MAE para los modelos incluidos en la lista.
map_dfr(.x = lista_pred_training, .f = metrics, truth = peso, estimate = .fitted, .id="modelo") %>% arrange(.estimate)

Luego el modelo que tiene menor MAE (7.3946489) y RMSE (9.8198125) sigue siendo el modelo propuesto 1 y el que le sigue es el modelo categóricas

En testing

#Modelo inicial
pred_modelo_ini_test = augment(modelo_inicial, newdata = dt_test)
metrics(data = pred_modelo_ini_test, truth = peso, estimate= .fitted)

#Modelo categóricas con variable con categoría __consumoSemanalSnacks__ redefinida
pred_modelo_interac_2_test = augment(modelo_interacion_2, newdata = dt_test_2)
metrics(data = pred_modelo_interac_2_test, truth = peso, estimate= .fitted)

#Modelo propuesto 1
pred_modelo_propuesto1_test = augment(modelo_propuesto_1, newdata = dt_test)
metrics(data = pred_modelo_propuesto1_test, truth = peso, estimate= .fitted)

#Modelo propuesto 2
pred_modelo_propuesto2_test = augment(modelo_propuesto_2, newdata = dt_test_propuesta_2)
metrics(data = pred_modelo_propuesto2_test, truth = peso, estimate= .fitted)

El MAE más bajo lo tiene el modelo propuesto 1 con 7.5561039 (pero con el modelo siguiente modelo categóricas solo tiene una diferencia de 0.0013), mientras que el RSME más bajo lo tiene el modelo categóricas con 10.1861909.

#Cálculo de R cuadrado ajustado (https://statologos.jaol.net/r-ajustado-al-cuadrado-en-r/)
adjusted_r2 <- function (modelo,r.squared,nro_obs) {
    return (1 - ((1-r.squared) * ( nro_obs -1) / (nro_obs - length(modelo$coefficients) -1)))
}

adjusted_r2(modelo_inicial, 0.3433755, nrow(dt_test))

## [1] 0.342064

adjusted_r2(modelo_interacion_2, 0.3467361, nrow(dt_test_2))

## [1] 0.344777

adjusted_r2(modelo_propuesto_1, 0.3428871, nrow(dt_test))

## [1] 0.3297493

adjusted_r2(modelo_propuesto_2, 0.3426815, nrow(dt_test_propuesta_2))

## [1] 0.3413686

En el caso del R² ajustado el valor de mayor variabilidad se obtiene con el modelo categóricas

Finalmente el mejor modelo para predecir el peso es el modelo categóricas , ya que fue el que mejores métricas devolvió utilizando el dataset de testing. Es decir, evaluando sobre datos que nunca habia visto.

Diagnóstico del modelo

Análisis de los supuestos del modelo lineal para el modelo inicial

plot(modelo_inicial)

Residuos versus valores predichos: No parece existir una estructura clara de datos, por lo que se podría decir que se cumple el supuesto de homocedasticidad.

Normal QQ plot: El extremo superior derecho no se ajusta a la distribución teórica (∼N(0,1)), por lo que no parece seguir una distribución normal.

Residual versus leverage: El leverage mide cuan influentes son los puntos en el modelo. En este gráfico algunas observaciones aparecen bastante separadas, indicando un alto leverage. El gráfico destaca las observaciones 6818, 2094 y 93. Significa que estas observaciones cambian mucho al modelo.

Diagnóstivo del modelo: El modelo no cumple con los supuestos del modelo lineal. Parecen existir 2 problemas: falta de normalidad y la presencia de observaciones con alto leverage.

Modelo robusto

Se incorporan al dataset original de train algunas observaciones adicionales que pueden incluir valores atípicos. Se quiere observar en particular la relación entre peso y altura:

ggplot(data = dt_train, aes(x = peso, y= altura)) +
  geom_point(alpha = 0.75) + 
  labs(title = "Peso y Altura (Modelo original)") + 
  theme(legend.position = 'none')

ggplot(data = dt_train_6, aes(x = peso, y= altura)) +
  geom_point(alpha = 0.75) + 
  labs(title = "Peso y Altura (Modelo con posibles outliers)") + 
  theme(legend.position = 'none')

Cambia la relación entre el peso y la altura. En el dataset con posibles outliers las personas (según peso y altura) se dividen en 2 grupos, uno más denso donde a mayor peso crece más rápidamente la altura que con el dataset original y otra agrupación adicional de puntos hacia la derecha con pesos más grandes que no existían en el dataset original.

Entrenamiento del modelo inicial con nuevos datos.

modelo_inicial_con_outliers = lm(formula = peso ~ altura + edad + genero + act_fisica_s + frec_alcohol_d, data= dt_train_6)
tidy_modelo_inicial_con_outliers  <- tidy(modelo_inicial_con_outliers, conf.int = TRUE)

Comparación de coeficientes (modelo inicial entrenado con nuevos datos versus con el set de datos original)

tidy_modelo_inicial_con_outliers

tidy_modelo_inicial

En cuanto a los coeficientes, si bien los valores estimados y sus errores standard cambian, no se producen grandes saltos. Además, siguen siendo estadísticamente significativos los mismos coeficientes en ambos modelos.

Comparación de métricas (modelo inicial entrenado con nuevos datos versus con el set de datos original)

pred_modelo_ini_con_outliers = augment(modelo_inicial_con_outliers, newdata = dt_train_6)
metrics(data = pred_modelo_ini_con_outliers, truth = peso, estimate= .fitted)

pred_modelo_ini_ori = augment(modelo_inicial, newdata = dt_train)
metrics(data = pred_modelo_ini_ori, truth = peso, estimate= .fitted)

El RMSE del modelo original es de 9.9101379 y el MAE es de 7.4552017, ambos menores que el del modelo con outliers.

glance(modelo_inicial_con_outliers)

glance(modelo_inicial)

El R² ajustado del modelo inicial logra explicar una variabilidad en un 35.4% mientras que el modelo con outliers solo logra explicar el 27.2%. Todas las métricas dan mejor para el modelo original

Entrenamiento con un modelo robusto.

library(robustbase)

modelo_inicial_robusto = lmrob(formula = peso ~ altura + edad + genero + act_fisica_s + frec_alcohol_d, data= dt_train_6)
tidy_modelo_inicial_robusto  <- tidy(modelo_inicial_robusto, conf.int = TRUE)

Comparación de coeficientes (Robusto y no Robusto)

tidy_modelo_inicial_con_outliers

tidy_modelo_inicial_robusto

Comparación de performance (Robusto y no Robusto)

pred_modelo_inicial_robusto = augment(modelo_inicial_robusto, newdata = dt_train_6)
metrics(data = pred_modelo_ini_con_outliers, truth = peso, estimate= .fitted)

metrics(data = pred_modelo_inicial_robusto, truth = peso, estimate= .fitted)

Con el ajuste robusto se pueden detectar outliers mirando aquellas observaciones cuyos pesos (robustos) asignados por el ajuste de lmrob sean muy chicos (pesos cero o muy cercanos a él).

summary(modelo_inicial_robusto)

## 
## Call:
## lmrob(formula = peso ~ altura + edad + genero + act_fisica_s + frec_alcohol_d, 
##     data = dt_train_6)
##  \--> method = "MM"
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -25.9082  -5.6143  -0.1514   6.0878 115.4825 
## 
## Coefficients:
##                  Estimate Std. Error t value Pr(>|t|)    
## (Intercept)     -63.71242    2.25261 -28.284  < 2e-16 ***
## altura            0.62359    0.01449  43.041  < 2e-16 ***
## edad              1.27269    0.08694  14.638  < 2e-16 ***
## generoMasculino   1.11503    0.25240   4.418 1.01e-05 ***
## act_fisica_s     -0.00211    0.04592  -0.046    0.963    
## frec_alcohol_d    0.07219    0.05579   1.294    0.196    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Robust residual standard error: 8.579 
## Multiple R-squared:  0.3869, Adjusted R-squared:  0.3865 
## Convergence in 12 IRWLS iterations
## 
## Robustness weights: 
##  62 observations c(93,248,779,1200,1389,1657,1702,2094,2232,2637,2882,2998,3002,3656,4084,4105,4631,5043,5174,5813,6114,6251,6369,6683,6818,6956,7025,7026,7027,7028,7029,7030,7031,7032,7033,7034,7035,7036,7037,7038,7039,7040,7041,7042,7043,7044,7045,7046,7047,7048,7049,7050,7051,7052,7053,7054,7055,7056,7057,7058,7059,7060)
##   are outliers with |weight| = 0 ( < 1.4e-05); 
##  572 weights are ~= 1. The remaining 6426 ones are summarized as
##      Min.   1st Qu.    Median      Mean   3rd Qu.      Max. 
## 0.0000908 0.8667000 0.9513000 0.8928000 0.9858000 0.9990000 
## Algorithmic parameters: 
##        tuning.chi                bb        tuning.psi        refine.tol 
##         1.548e+00         5.000e-01         4.685e+00         1.000e-07 
##           rel.tol         scale.tol         solve.tol       eps.outlier 
##         1.000e-07         1.000e-10         1.000e-07         1.416e-05 
##             eps.x warn.limit.reject warn.limit.meanrw 
##         3.638e-10         5.000e-01         5.000e-01 
##      nResample         max.it       best.r.s       k.fast.s          k.max 
##            500             50              2              1            200 
##    maxit.scale      trace.lev            mts     compute.rd fast.s.large.n 
##            200              0           1000              0           2000 
##                   psi           subsampling                   cov 
##            "bisquare"         "nonsingular"         ".vcov.avar1" 
## compute.outlier.stats 
##                  "SM" 
## seed : int(0)

En este caso se detectaron 62 observaciones como outliers.

Ahora bien, si comparamos las métricas de ambos modelos, el RMSE es mayor en 0.0005 y el MAE es mayor en 0.1115 en el modelo con outliers respecto del modelo robusto. Es decir que las medidas de performance dan casi similares, es solo infimamente mejor en el modelo robusto. Entonces, en este caso se puede concluir que la presencia de outliers no afectó demasiado al modelo original (que es no robusto).

TP 1: Regresión Lineal

Claudia Roxana Aguirre

29 de Octubre, 2021

Preparación de datos

Análisis exploratorio

Estructura y variables

Correlación entre variables numéricas

Apertura por género

Hambre mensual, consumo semanal de verdura y consumo semanal de comida grasa

Modelo Inicial

Coeficientes estimados

Interpretación

Significatividad

Modelo categóricas

Coeficientes estimados

Interpretación

consumoSemanalSnacks

genero.edad

Variable consumoSemanalSnacks

Evaluación de significatividad conjunta

Redefinición de las categorías:

Modelos propios y evaluación

Creación de 2 modelos múltiples adicionales

Modelo 1 propuesto

Modelo 2 propuesto

Evaluación de la performance de los modelos

En training

En testing

Diagnóstico del modelo

Modelo robusto

Entrenamiento del modelo inicial con nuevos datos.

Comparación de coeficientes (modelo inicial entrenado con nuevos datos versus con el set de datos original)

Comparación de métricas (modelo inicial entrenado con nuevos datos versus con el set de datos original)

Entrenamiento con un modelo robusto.

Comparación de coeficientes (Robusto y no Robusto)

Comparación de performance (Robusto y no Robusto)